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Abstract A new law of motion for supernova remnant 
(SNR) which introduces the quantity of swept matter 
in the thin layer approximation is introduced. This new 
law of motion is tested on 10 years observations of SN 
1993J . The introduction of an exponential gradient in 
the surrounding medium allows to model an aspherical 
expansion. A weakly asymmetric SNR, SN 1006 , and 
a strongly asymmetric SNR, SN 1987A , are modeled. 
In the case of SN 1987A the three observed rings are 
simulated. 

Keywords supernovae: general supernovae: individ- 
ual (SN 1993J ) supernovae: individual (SN 1006 ) 
supernovae: individual (SN 1987A ) ISM : supernova 
remnants 



1 Introduction 

The theoretical study of supernova remnant (SNR) has 
been focalized on an expression for the law of mo- 
tion. As an example the Sedov- Taylor expansion pre- 
dicts R cx t°-^, see see 



than the Sedov- Taylor solution, more details for the 



Taylor (1950); Sedov (1959); 



McCray and Layzer (1987) and the thin layer approx- 



imation in the presence of a constant density medium 
predicts R oc see Dyson, J. E. and Wilhams, D. A. 



(1997); Dyson (1983); [Canto, Raga, and Adame (2006). 



The very-long-baseline interferometry (VLBI) observa- 
tions of SN 1993J (wavelengths of 3.6, 6, and 18 cm) 
show that R cx t^'^^ over a 10 year period, see Mar- 



caide et al. (2009). This observational fact does not 



agree with the current models because the radius of SN 
1993J grows slower than the free expansion and faster 
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spherical case can be found in Zaninetti (2011). The 



SNRs can also be classified at the light of the observed 
symmetry. A first example is SN 1993J which presented 
a circular symmetry for 4000 days, see [Marcaide et al\ 
( 2009 1 . An example of weak departure from the circu- 



lar symmetry is SN 1006 in which a ratio of 1.2 be- 
tween maximum and minimum radius has been mea- 



sured, see Reynolds and Gilmore (1986). An example 



of axial symmetry is SN 1987A in which three rings are 
symmetric in respect to a line which connect the cen- 



ters, see Tziamtzis et al. (2011 ). The models cited leave 



some questions unanswered or only partially answered: 

• Is it possible to deduce an equation of motion for an 
expanding shell assuming that only a fraction of the 
mass enclosed in the advancing sphere is absorbed in 
the thin layer? 

• Is it possible to model the complex three-dimensional 
(3D) behavior of the velocity field of the expanding 
nebula introducing an exponential law for the density 
? 

• Is it possible to make an evaluation of the reliabil- 
ity of the numerical results on radius and velocity 
compared to the observed values? 

• Can we reproduce complicate features such as equa- 
torial ring -(- two outer rings in SN 1987A which are 
classified as a "mystery" ? 

• Is it possible to build cuts of the model intensity 
which can be compared with existing observations? 

In order to answer these questions. Section [2] de- 
scribes two observed morphologies of SNRs, Section [373] 
reports a new classical law of motion which introduces 
the concept of non cubic dependence (NCD) for the 
mass included in the advancing shell. Section |4] intro- 
duces an exponential behavior in the number of parti- 
cles which models the aspherical expansion, Section [5] 
applies the law of motion to SN 1987A and SN 1006 
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introducing the quality of the simulation, Section |6] re- 
views the existing situation with the radiative transport 
equation and Section[7]contains detailed information on 
how to build an image of the two astrophysical objects 
here considered. 



2 Astrophysical Objects 

2.1 A strongly asymmetric SNR , SN 1987A 

The SN 1987A exploded in the Large Magellanic Cloud 
in 1987. The distance of this SN is « 50 kpc (163050 ly) 
and a detailed analysis of the distance , D , gives D = 



3 Classical law of motion 



3.1 Thermal or non thermal emission? 



The synchrotron emission in SNRs is detected from 
lO^Hz of radio-astronomy to lO^^Hz of gammaastrono- 
my which means 11 decades in frequency. At the same 
time some particular effects such as absorption, tran- 
sition from optically thick to optically thin medium, 
line emission, and energy decay of radioactive isotopes 
( ^^Ni, ^^Co) can produce a change in the concavity of 
the flux versus frequency relationship, see the discus- 
sion about Cassiopea A in Section 3.3 of |Eriksen et al.\ 
{ 2009 1 . A comparison between non-thermal and ther- 



51.4 kpc Panagia (2005) and D — 50.18 kpc Mitchell mal emission (luminosity and surface brightness distri- 



ct al. (2002). In the numerical codes we will assume 
D = 50kpc . The observed image is complex and we 



will follow the nomenclature of Racusin et al. (2009) 



which distinguish between torus only , torus +2 lobes 
and torus -I- 4 lobes. In particular we concentrate on 
the torus which is characterized by a distance from the 
center of the tube and the radius of the tube. In Table 



2 of Racusin et al. { 2009 ) is reported the relationship 



between distance of the torus in arcsec and time since 
the explosion in days. 

2.2 A weakly asymmetric SNR, SN 1006 

This SN started to be visible in 1006 AD and actually 
has a diameter of 12.7 pc, see Strom ([1988^) . More 
precisely, on referring to the radio-map of SN 1006 at 



1370 MHz by Reynolds and Gilmore (1986), it can be 



observed that the radius is greatest in the north-east 
direction. From the radio-map previously mentioned 
we can extract the following observed radii, R — 6.8 pc 
in the polar direction and R = 5.89 pc in the equato- 
rial direction. Information on the thickness of emitting 



layers is contained in Bamba et al. (20031 where the 



Chandra observations (i.e., synchrotron X-rays) from 
SN 1006 were analyzed. The observations found that 
sources of non-thermal radiation are likely to be thin 
sheets with a thickness of about 0.04 pc upstream and 
0.2 pc downstream of the surface of maximum emission, 
which coincide with the locations of Balmer-line optical 
ilEllison et a/.|(|1994|) . The high resolution 



emission , see 



XMM-Newton Reflection Grating Spectrometer (RGS) 
spectrum of SN 1006 gives two solutions for the O VII 
triplet. One gives a shell velocity of 6500 km/ s and 
the second one a shell velocity of 9500 km/s, when a 



distance of 3.4 kpc is adopted ,see Vink (2005) 



bution) can be found in Petruk and Beshlei (2007), 
where it is possible to find some observational tests 
which allow estimation of parameters characterizing the 
cosmic ray injection on supernova remnant shocks. At 
the same time, a technique to isolate the synchrotron 
from the thermal emission is widely used, as an example 



see X-limb of SN1006 Katsuda et al. (2010). 



3.2 The basic assumptions 

The observational dichotomy between thermal and non 
thermal emission influences the use of temperature in 
the theoretical models. As an example the Sedov solu- 
tion is 



R(t) « 0.313 



E51 t 



51 1-1 

no 



pc 



(1) 



where ti is the time expressed in years, -B51, the energy 
in 10^^ erg and tlq is the number density expressed in 
particles cm~'^ (density p = nom, where m = 1.4mH), 
see|Sedov|(|1959|);|McCray and Layzer|(|1987|);|Zaninetti| 



(2011 1. The spectrum of the emitted radiation depends 
on the temperature behind the shock front , see for 



example formula 9.14 in McKee (1987) 



1 = K 

16 A: ^ 



(2) 



where p, is the mean mass per particle, k the Boltz- 
mann constant and vs the shock velocity expressed in 
cm sec~^. On identifying Vs with the velocity derived 
from Eq. [TJthe Sedov solution has the following behav- 
ior in temperature 



no'^/^ti^'^ 



(3) 



when p, is the mass of the hydrogen. The Sedov ap- 
proach cannot be applied to the aspherical case because 
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we do no know the complex volume occupied during the 
expansion. A numerical approach developed to follow 
the evolution of the superbubbles first compute the vol- 
ume occupied during the expansion and then the pres- 
sure , see Zaninetti (2004). Conversely the thin layer 



see 



Dyson, J. E. and Williams, D. A.| 



approximation 

([l997j); ,Dyson| ( |1983[ ), takes in consideration only the 
swept mass and the velocity is deduced applying the 
momentum conservation . Due to the absence of the 
temperature and pressure the thin layer approximation 
can be classified as a non thermal model. 



is: 



where 
N = 

(i?oi+^^'"+(3+p)x 



D = Ro 



l+3p- 



3p- 



tp 



-SVoRo''''' to 



VoRo 



3p- 



top. 



(10) 



3.3 The incomplete thin layer approximation 

The thin layer approximation with non cubic depen- 
dence (NCD) , p , in classical physics assumes that 
only a fraction of the total mass enclosed in the volume 
of the expansion accumulates in a thin shell just after 
the shock front. The global mass between and i?o is 
^ttpRq where p is the density of the ambient medium. 
The swept mass included in the thin layer which char- 
acterizes the expansion is 



Mo^i-npRl)-. . 

The mass swept between and R is 



M^i-npR^ 



(5) 



The conservation of radial momentum requires that, 
after the initial radius Rq, 



MV = MoVo 



(6) 



where R and V are the radius and velocity of the ad- 
vancing shock. In classical physics, the velocity as a 
function of radius is: 



and introducing f3o — — and (5 — — we obtain 



/3 = /3o(f)i . 

The law of motion is: 



(7) 



(8) 



i?o'+'^" +{3+p)voRo^P' {t-to)p 
where t is time and Iq is the initial time. 



(9) 



Equation ^ can also be solved with a similar solution 
of type R = Kit — to)", k being a constant, and the 
classical result is: 



m = 



' ('3+p)VoRo''^P (t-to)^ 



(11) 



The similar solution for the velocity is 



V{t) = 



(4) (3+p) p^+p Vo'-'+p Ro^'^" {t - to) 



(12) 



The similar formula ( 11 ) for the radius can be com- 



pared with the observed radius radius-time relationship 
of the Supernova reported as 



R{t) = robst° 



(13) 



where the two parameters rof,s and a^bs are found from 
the numerical analysis of the observational data. In this 
case the velocity is 



Vit) = robsaob.t<'''"''-^^ . 



(14) 



The comparison between theory and astronomical ob- 
servations allows to deduce p , the NCD parameter , as 



3a, 



P = 



ohs 



1 - Olohs 



(15) 



where Uobs characterizes the radius-time relationship 
in SNRs as given by Eq. (13). The Supernova SN 



1993 J represents a test for this theoretical model of the 
expansion. A careful analysis of the radius-time rela- 
tionship for SN 1993J shows that R{t) sa O.OlStO-^^^^ 
when the time is expressed in yr . A first possible 
system of units which allows to make a comparison 
with the observations is represented by pc for the length 
and yr for the time. The theoretical solution as given 
by equation ^ can be found through the Levenberg- 
Marquardt method ( subroutine MRQMIN in Press 
et al. ( 1992 )) and Fig. [T] reports a numerical example. k 
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Time (yr) 

Fig. 1 Theoretical radius as obtained by the solution of the equation ([9} (full line), data and merit function , , as in 
Table [T] The astronomical data of SN 1993 J are represented through empty stars. 



The quality of the fits is measured by the merit func- where 



tion X 



E 



{Rth — RobsY 



(16) 



where Rth, Robs and aobs are the theoretical radius, 
the observed radius and the observed uncertainty re- 
spectively. 

The conservation of classical momentum here adop- 
ted does not take into account the momentum carried 
away by photons or in other words the radiative losses 
are included in the NCD exponent. 



4 Asymmetrical law of motion with NCD 

Given the Cartesian coordinate system (x, y, z) , the 
plane z — will be called equatorial plane and in polar 
coordinates z — Rsm{9), where 9 is the polar angle 
and R the distance from the origin . The presence of a 
non homogeneous medium in which the expansion takes 
place is here modeled assuming an exponential behavior 
for the number of particles of the type 



n{z) — no exp - 



no exp 



R X sin(6') 



(17) 



where R is the radius of the shell, uq is the number of 
particles at R = Rq and h the scale. The 3D expansion 
will be characterized by the following properties 

• Dependence of the momentary radius of the shell on 
the polar angle that has a range [—90° O -1-90°]. 

• Independence of the momentary radius of the shell 
from (f> , the azimuthal angle in the x-y plane, that 
has a range [0° O 360°]. 

The mass swept, M, along the solid angle A fl, between 
and R is 

M{R) = ^mnnoImiR) + ^Tri^gnom^/ , (18) 



Irn{R)^ I r 2 exp - ' "7'"' dr 



Ro 



r sin(6') 



(19) 



where Rq is the initial radius and mn the mass of the 
hydrogen. The integral is 



ImiR) = 



h (2h^ + 2Rohsm (9) + Rq^ (sin {6)f^ e" 



(sin(0))^ 



h[2h^ + 2Rh sin {6) + R^ (sin {Oyy e 



2\ H3in(9) 



(20) 



The conservation of the momentum gives 
{M{R))^R^{M{Ro))"^Ro , (21) 

where R is the velocity at R and Rq is the initial velocity 
at R — Rq where p is the NCD parameter. This means 
that only a fraction of the total mass enclosed in the 
volume of the expansion accumulates in a thin shell 
just after the shock front. According to the previous 
expression the velocity is 



ViR) = Ro'P 'Vo{i-2h^e 



„ R sin 



-2h^e ^i?sin(6l) 



-he —R^{siii{9)f + 2h^e~ 

+2h^e s — Rg sm{e) + he s — Rg^ {fim{e)f 



(22) 



In this differential equation of the first order in R the 
variable can be separated and the integration term by 
term gives 



(M(r))i/f dr = {M{Ro)y^PRo x (t - to) 



(23) 



-Ro 
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where t is the time and to the time at Rq. The result- 
ing non hnear equation Tnl expressed in astrophysical 
units is 



NL 



{{-2h 



n(e) 



pc ^ 

n(6) 



' pc 



sin(6l) 



d„(8) 



'Zhp^ e 



-he '^p" rpc {sm{6)) 



"p- i?o,pcSin(6l) 



+/ie i?o,pc'(sin(0))2 
+ii;o,pc'(sin(0))3)(sin(0))-3)p-'dr 



O.pc'^^ Ro,kms{tl — to.l) — (24) 



where ii and to,i are t and to expressed in yr units, rpc 
and Ro.pc are i? and i?o expressed in pc, Ro^kms is Ro 
expressed in is expressed in radians and hpc is 
the the scale , h , expressed in pc. It is not possible to 
find Rpc analytically and a numerical method should 
be implemented. In our case in order to find the root 
oiJ^NL, the FORTRAN SUBROUTINE ZRIDDR from 



Press et at. ( 1992 1 has been used 



5 Applications of the law of motion 

From a practical point of view, e , the percentage of 
reliability of our code is introduced, 



(1 



K^pc 



obs 



100, 



(25) 



where Rpc,ohs is the radius as given by the astronomical 
observations in parsec , and i?pc,num the radius obtained 
from our simulation in parsec. 



Table 1 Numerical value of the parameters of the fits 
for SN 1993 J and x^- N represents the number of free pa- 
rameters. PL stands for power law and NLR for non linear 
radius. 



5.1 Results on the strongly asymmetric SN 198 7 A 

The first target is to simulate the torus only of SN 
1987A (our equatorial plane) and this operation allows 
to calibrate our code ; Table [2] reports the input data. 
Figs. [2]and[3]the behavior of radius and velocity respec- 
tively as function of the time. Due to the complexity of 
the 3D structure we report the reliability of the simula- 
tion only in the equatorial plane (torus only) which is 
e = 91.4% for the radii and e = 38.3% for the velocity. 

After this calibration on the equatorial plane we con- 
tinue identifying the lobes of SN 1987A as bipolar SNR 
as seen from a given point of view. The complex 3D be- 
havior of the advancing SNR is reported in Fig. [4] and 
Fig. [5] reports the asymmetric expansion in a section 
crossing the center. In order to better visualize the 
asymmetries Figs. |6] and [7] report the radius and the 
velocity respectively as a function of the position angle 
9. The combined effect of spatial asymmetry and field 
of velocity are reported in Fig. [Sj 

An explanation for the lack of velocity in the equa- 
torial region of SN 1987A can be drawn from a careful 
analysis of Fig. [9] which displays the swept mass as a 
function of the latitude. The ratio between maximum 
swept mass in the polar direction and minimum swept 
mass in the equatorial plane is « 5.5 . On applying the 
momentum conservation the velocity in the equatorial 
plane is 5.5 times smaller in respect to the polar direc- 
tion where the velocity is maximum. This numerical 
evaluation gives a simple explanation for the asymme- 
try of SN 1993J : a smaller mass being swept means a 
greater velocity of the advancing radius of the nebula. 

5.2 Results on the weakly asymmetric SN 1006 

The input data of the simulation are reported in Table|3] 
; the reliability is e = 89% for the radii in the equatorial 
direction and e = 89% for the velocity in the equatorial 
direction. 

The weakly asymmetric 3D shape of SN 1006 is re- 
ported in Fig. [To] 



keyword 


N 


values 




PL 


2 


aobs = 0.82; robs = 0.015 pc 


6364 


NLR 


4 


p^ 18.8; 


623.5 



ro = 2.25 10-^ pc; 
to = 4.98 0-*^ yr; 



100000^ 



Table 2 Numerical value of the parameters of the simu- 
lation for SN 1987A 



Ro,pc 

RQ,kms 
P 

*0,1 

h 



0.014 
26000 
4 

0.0218 
21.86 
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Fig. 2 Radius as a function of time for an exponen- 
tially varying medium (full line) when 9 = 0.001 (equatorial 
plane) and astronomical data of torus only as extracted from 
Table 2 in |Racusin et aZ.] ( |2009[ ). Physical parameters as in 
Tabled 




Time 



Fig. 3 Velocity as a function of the time for an exponen- 
tially varying medium (full line) when 9 = 0.001 (equatorial 
plane) and astronomical velocities as deduced from Table 2 
in Racusin et al. (20091. Physical parameters as in Table [2] 




Fig. 4 Continuous three-dimensional surface of SN 1987A 
: the three Eulerian angles characterizing the point of view 
are $=90 °, 9=90 ° and ^'=90 °. Physical parameters as 
in Table [21 




X [pc] 



Fig. 5 Section of SN 1987A in the x-z plane. The hor- 
izontal and vertical axis are in pc. Physical parameters as 
in Table H 




d (degree) 



Fig. 6 Radius in pc of SN 1987A as a function of the 
position angle in degrees. Physical parameters as in Table 

m 




(degree) 



Fig. 7 Velocity in ^ of SN 1987A as a function of the 
position angle in degrees. Physical parameters as in Table 

m 
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25003.510 km/s«c 




The velocity as function of the position angle is plot- 



ted in Fig. 11 and a comparison should be done with 



Fig. 4 in Katsuda et al. ( 2009 1 where the proper motion 



as function of the azimuth angle was reported. 

Our model for SN 1006 predicts a minimum velocity 
in the equatorial plane of 5785 km/ s and a maximum 
velocity of 7229 km/ s in the polar direction. A recent 
observation of SN 1006 quotes a minimum velocity of 
5500 km/s and a maximum velocity of 14500 km/s 
assuming a distance of = 3.4 kpc The swept mass in 
the thin layer versus the position angle is displayed in 
Fig. Ill 




5400.791 km/sac 

Fig. 8 Map of the expansion velocity in ^ relative to 
the simulation of SN 1987A when 300000 random points 
are selected on the surface. Physical parameters as in Table 

m 




Fig. 9 Swept mass of SN 1987A as a function of the posi- 
tion angle in degrees in arbitrary units. Physical parameters 
as in Tabled 



Table 3 Numerical value of the parameters of the simu- 
lation for SN 1006 . 



Ro,pc 
Ro,kr, 
P 

^0.1 

h 



0.09 
8500 

30 

10 
1000 



5.3 Comparison with the stellar wind 

A comparison can be done with the expansion speed of 
the outer H2 shell of 77-Carinae which has been fitted 
with the following latitude dependent velocity 



vi (1 + e2-^ cos(2e)) 



(26) 



where the parameter A controls the shape of the Ho- 
munculus, Q is the polar angle; vi and V2 are the veloc- 



ities in the polar and equatorial direction, see Gonzalez 



et al. (2010). Fig. 13 the data of our simulation as well 



the wind-type profile of velocity. 



6 Radiative transfer equation 

The transfer equation in the presence of emission only 



ing, R. M. (1988), is 



see for example Rybicki and Lightman ( 1991 ) or Hjellm- 



ds 



(27) 



where is the specific intensity , s is the line of sight 
, ji, the emission coefficient, ki, a mass absorption coef- 
ficient, C, the mass density at position s and the index 
V denotes the interested frequency of emission. The 



solution to equation (27) is 



/,(r.) = Z^(l_e--^(^)) , (28) 

where r,^ is the optical depth at frequency u 

dr^ = k^Cds . (29) 

We now continue analyzing the case of an optically thin 
layer in which t^, is very small ( or k^ very small ) and 
the density C, is substituted with our number density 
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Fig. 10 Continuous tliree-dimensional surface of tiie SN 
1006 ; the three Eulerian angles characterizing the point of 
view are $=90 °, 0=90 ° and ^=90 °. Physical parameters 
as in Tabled 




(degree) 



Fig. 11 Velocity in ^ of SN 1006 as a function of the 
position angle in degrees. Physical parameters as in Table 



H 



C(s) of particles. One case is taken into account : the 
emissivity is proportional to the number density. 



= KC{s) 



(30) 



where X is a constant function. This can be the case 
of synchrotron radiation in presence of a isotropic dis- 
tribution of electrons with a power law distribution in 
energy, N{E), 

N{E)dE = KsE-'^f , (31) 
where Ks is a constant. In this case the emissivity is 



0.933 X 10~2^a(7/)if^i7^ 
.6.26 X lO^^Ii.:! erg 



21+1 



scm^ H zracP 



(32) 



where v is the frequency and 0^(7/) is a slowly varying 
function of 7/ which is of the order of unity and is given 

by 



g(-y,-3)/2 7/ + 7/3 . 37/ - ^ p^37/ + 
7/ + 1 ^ 12 12 



a (7/) 
7. 



(33) 



for 7/ > 5, see formula (1.175 ) in Lang (1999). The 



source of synchrotron luminosity is assumed here to be 
the flux of kinetic energy, Lm , 



L„, = ^pAV , (34) 
where A is the considered area, see formula (A28) in de 



Young (2002). In our case A = AttR^, which means 
1 



,2ta3 



(35) 



where R is the instantaneous radius of the SNR and 
p is the density in the advancing layer in which the 




100 200 300 

S (degree) 




20 40 60 80 

8 (degree) 



Fig. 12 Swept mass of SN 1006 as a function of the posi- 
tion angle in degrees in arbitrary units. Physical parameters 
as m Tabled 



Fig. 13 Expansion velocity versus latitude of SN 1987A 
( open stars) and fit as given by formula (261 (full line). 
The fitting parameters are : A=1.9, vi = 2503.6 fcm/s and 
V2 = 5401 km/ s. 
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synchrotron emission takes place. The total observed 
luminosity can be expressed as 



^tot — ^Lra 



(36) 



where e is a constant of conversion from the mechani- 
cal luminosity to the total observed luminosity in syn- 
chrotron emission. The fraction of the total luminosity 
deposited in a given band fb is 



^ , min '^b,max 

y . -7/ + 1 _ y -7/ + 1 



(37) 



where Vb.min and Vb,max are the minimum and maxi- 
mum frequency of the given band. 



7 Image 

An simulated image of an astrophysical object is com- 
posed by combining the intensities which character- 
ize different points. For an optically thin medium the 
transfer equation provides the emissivity to be multi- 
plied with the distance on the line of sight , /. This 
length in astrophysical diffuse objects depends on the 
orientation of the observer. A thermal and a non ther- 
mal model are reviewed in the spherical case. In the 
the aspherical case a non thermal model is presented. 

7.1 Spherical Image 

A first thermal model for the image is characterized by 
a constant temperature in the internal region of the ad- 
vancing sphere. We therefore assume that the number 
density C is constant in a sphere of radius a and then 
falls to 0. The length of sight , when the observer is sit- 
uated at the infinity of the x-axis , is the locus parallel 
to the a;-axis which crosses the position y in a Carte- 
sian X — y plane and terminates at the external circle 



of radius a, see Zaninetti (2009). The locus length is 



y2) ;0 < y < a 



(38) 



The number density is constant in the sphere of 
radius a and therefore the intensity of radiation is 



/oa = X 2 X (Va 



■,0<y <a 



(39) 



The comparison of observed data of SN 1993J and the 



theoretical thermal intensity is reported in Fig. 14 



A second non thermal model for the image is charac- 
terized by emission in a thin layer around the advancing 
sphere. We therefore assume that the number density 
C is constant and in particular rises from at r = a to 
a maximum value , remains constant up to r = 6 



and then falls again to 0. The length of sight , when 
the observer is situated at the infinity of the x-axis , is 
the locus parallel to the x-axis which crosses the posi- 
tion y in a Cartesian x — y plane and terminates at the 
external circle of radius 6, see 
locus length is 



Zaninetti (2009). The 



^oa = 2 X ( - 2/2 _ _ y2) ;0<y<a 

■,a<y <b . (40) 



U = 2 X i^b^ - y2) 



The number density Cm is constant between two 
spheres of radius a and b and therefore the intensity 
of radiation is 



/oa = C„ X 2 X (v/fe2 - y2 _ _ y2) ;0<y<a 

Iab^CmX2x ( - y2) ;a<y<b .(41) 

The comparison of observed data of SN 1993J and the 
theoretical non thermal intensity is displayed in Fig. 

The main result of this Section is that the intensity 
of the thermal model which has the maximum of the 
intensity at the center of SNR does not match with the 
observed profiles. The observed profiles in the intensity 
have the maximum value at the rim as predicted by the 
non thermal model. 

7.2 Aspherical Image 

The numerical algorithm which allows us to build a 
complex image is now outlined. 

• An empty (value^O) memory grid Ai{i,j,k) which 
contains NDIM^ pixels is considered 

• We first generate an internal 3D surface by rotat- 
ing the ideal image 180° around the polar direction 
and a second external surface at a fixed distance AR 
from the first surface. As an example, we fixed AR 



0.03i?„ 



where Rmax is the maximum radius of 



expansion. The points on the memory grid which lie 
between the internal and external surfaces are mem- 
orized on A4(i,j,k) with a variable integer number 
according to formula ( |35[ ) and density p proportional 
to the swept mass, see Fig. |9] 
• Each point of A^(z, j, k) has spatial coordinates x, y, z 
which can be represented by the following 1x3 ma- 
trix, A, 



A = 



(42) 



The orientation of the object is characterized by the 
Euler angles ($, G, ^') and therefore by a total 3x3 



10 



JH^^ 
















\ 





Radius (mas) 



Fig. 14 Cut of the thermal intensity / of the rim model ( equation (39l) through the center (dotted line) of SN 1993 J and 
real data (empty stars) when 6 = 0.07pc. The observed data as day 1889 after the explosion have been extracted by the 
author from Fig. 3 of|Marcaide et a/.| ([20091. 























% 





Radius (mas) 



Fig. 15 Cut of the non thermal intensity / of the rim model ( equation (41l) through the center (dotted line) of SN 
1993Jand real data (empty stars). The parameters are a = 0.0035 pc and b = 0.07pc. Observed data as day 1889 after the 
explosion have been extracted by the author from Fig. 3 of Marcaide et al. (2009) . 
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rotation matrix, E, see Goldstein, Poole, and Safko 



( 2002 1 . The matrix point is represented by the fol- 



lowing 1x3 matrix, B, 



B = E- A 



(43) 



• The intensity map is obtained by summing the points 
of the rotated images along a particular direction. 

• The effect of the insertion of a threshold intensity , 
Itr, given by the observational techniques , is now an- 
alyzed. The threshold intensity can be parametrized 
to Imaxj the maximum value of intensity character- 
izing the map. 

7.3 The image of the strongly asymmetric SN 1987A 

An ideal image of SN 1987A having the polar axis 
aligned with the z-direction which means polar axis 
along the z-direction, is shown in Fig. |16| A model 
for a realistically rotated SN 1987A is shown in Fig. 

HH 

The three rings of SN 1987A are now simulated in 
our Fig. [it] and a comparison should be done with Fig. 
1 in [Tziamtzis et al.\ ( |2011[ ) in which the HST/ACS 
image (filter F658N) of the triple ring system of SN 
1987A was reported. 




I 



I 



-0.5 0.5 

orcsec from center 



Fig. 16 Map of the theoretical intensity of SN 1987A in 
the presence of an exponentially varying medium. Physical 
parameters as in Table |2] The three Euler angles character- 
izing the orientation are $=180°, 9=90 ° and *=0 °. This 
combination of Euler angles corresponds to the rotated im- 
age with the polar axis along the z-axis. 



7.4 The image of the weakly asymmetric SN 1006 
The image of SN 1006 is visible in different bands such 



as radio, see 


Reynolds and Gilmore 


( 


19931; 


Reynoso 


(2007 


), optical , see 


Long 


(2007 


and X-ray , see 


Dyer, 


Reynolds, and Borkowski 


2004) 


; Katsuda et al. 


2010 





The 2D map in intensity of SN 1006 is visible in Fig. 
[I8l The intensity along the equatorial and polar direc- 
tion of our image is reported in Fig. [T9j a comparison 
should be done with Fig. 4 in [Dyer, Reynolds, and| 
Borkowski (2004). The projected flux as a function of 



the position angle is another interesting quantity to plot 



, see Fig. 20 and a comparison should be done with Fig. 



5 top right in Rothenflug et al. (2004) 



After the previous graphs is more simple to present 
a characteristic feature as the "jet appearance" visible 
in some maps , see our Fig. |21[ a comparison should 
be done with the X-map at 6.33-6.53 kev band visible 



in Fig. 3b by Yamaguchi et al. ( 2008 1 



8 Conclusions 

Law of motion 

We have deduced a new law of motion in spherical 
symmetry ( constant density) for an advancing shell 



-0.5 0.5 

pc from center 



I 



i 



Fig. 17 Model map of SN 1987A rotated in accor- 
dance with the observations, for an exponentially varying 
medium. Physical parameters as in Table [2] The three Eu- 
ler angles characterizing the orientation of the observer are 
'I>=105°, 9=55 ° and *=-165 °. This combination of Eu- 
ler angles corresponds to the observed image. In this map 
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-10 -5 „ <0 X 5 10 

R Cpc) 

Fig. 19 Two cut along perpendicular lines of / for the non rotated image of SNR SN 1006 . 




100 200 300 

Position angle (degrees) 



Fig. 20 Intensity as function of the position angle in degrees for SNR SN 1006 . 




I 



I 



orcsec from center 



Fig. 18 Non rotated map of SN 1006 for an exponentially 
varying medium. Physical parameters as in Table ID 




I 



-5 5 

pc from center 



Fig. 21 Model map of SN 1006 rotated in accor- 
dance with the X observations, for an exponentially varying 
medium. Physical parameters as in Table H] The three Eu- 
ler angles characterizing the orientation of the observer are 
<E>=90 °, e=-55 ° and VE'=-180 °. This combination of Eu- 
ler angles corresponds to the observed image. In this map 
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assuming that only a fraction of the mass which re- 
sides in the surrounding medium is accumulated in the 
advancing layer, see equation (|9|. The presence of an 
exponential law for the density transforms the spherical 
symmetry in axial symmetry and allows the appearance 
of the so called "bipolar motion", see the nonlinear as- 



trophysical equation ( 24 1 . 
Images 

The emissivity in the advancing layer is assumed to 
be proportional to the flux of kinetic energy, see equa- 



tion (34 1 where the density is assumed to be propor- 
tional to the swept material. This assumption allows 
to simulate particular effects such as the triple ring sys- 
tem of SN 1987A , see Fig. [17] Another curious effect 
is the "jet appearance" visible in the weakly symmet- 

The jet/counter jet effect 



ric SN 1006 , see Fig. 21 



plays a relevant role in the actual research , see dis- 



cussion in Section 5.2 in Fesen et al. (20061 where the 



jet appearance is tentatively explained by the neutrino 
heating , seelWalder et al.\\2mb\ or by the MHD jet , 



see 



Takiwaki et al. (2004). Here conversely we explain 



the appearance of the jet by the addition of three effects 



• An asymmetric law of expansion due to a gradient 
in density in respect to the equatorial plane which 
produces an asymmetry in velocity . 

• The direct conversion of the flux of kinetic energy 
into radiation. 

• The image of the SNR as the composition of integrals 
along the line of sight. 

According to the previous three ingredients the neu- 
trino heating mechanism is not necessary. The Magneto 
Hydrodynamic (MHD) approach is supposed to act in 
a hidden way on scales smaller than than the thickness 
of the advancing layer in order to accelerate the elec- 
trons to relativistic energies. A careful calibration of 
the various involved parameters can be done when cuts 
in intensity are available. 
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